Phase separation frustrated by the long range Coulomb interaction I: Theory 
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We analyze the combined effect of the long range Coulomb (LRC) interaction and of surface energy 
on first order density-driven phase transitions between two phases in the presence of a compensating 
rigid background. In the coexistence region we study mixed states formed by regions of one phase 
surrounded by the other in the case in which the scale of the inhomogeneities is much larger than 
the interparticle distance. Two geometries are studied in detail: spherical drops of one phase into 
the other and a layered structure of one phase alternating with the other. In the latter case we 
find the optimum density profile in an approximation in which the free energy is a functional of 
the local density (LDA). It is shown that an approximation in which the density is assumed to be 
uniform (UDA) within each phase region gives results very similar to those of the more involved 
LDA approach. Within the UDA we derive the general equations for the chemical potential and the 
pressures of each phase which generalize the Maxwell construction to this situation. The equations 
are valid for a rather arbitrary geometry. We find that the transition to the mixed state is quite 
abrupt i.e. inhomogeneities of the first phase appear with a finite value of the radius and of the 
phase volume fraction. The maximum size of the inhomogeneities is found to be on the scale of a few 
electric field screening lengths. Contrary to the ordinary Maxwell construction, the inverse specific 
volume of each phase depends here on the global density in the coexistence region and can decrease 
as the global density increases. The range of densities in which coexistence is observed shrinks as the 
LRC interaction increases until it reduces to a singular point. We argue that close to this singular 
point the system undergoes a lattice instability as long as the inverse lattice compressibility is finite. 

Pacs Numbers: 64.75. Ig 71.10.Hf 71.10.Ca 



(N 
> 
(N 

a\ 
o 
o 

o 
o 



o 
o 



X 
J3 



I. INTRODUCTION 

The complex phase-diagrams of hole doped cuprates 
and manganites have rekiridled the study of mixed states 
in modeling these systemsEij Indeed strongly correlated 
systems with narrow bandwidth and short range interac- 
tions show a generic tendency to phase separate into hole- 
rich and hole-poor regions. When long range Coulomb 
(LRC) forces are taken into account this instability with 
macroscopic separation is frustrated due to the electro- 
static energy cost and this cariJead to charge inhomoge- 
neous states of various naturepu where domains of var- 
ious forms of one phase {B) are embedded in the other 
phase [A). 

In the inhomogeneous state the charge is segregated 
locally over some characteristic distance but the overall 
density (averaged over much larger distances) is a fixed 
constant in order to guarantee large scale neutrality and 
avoid the large Coulomb cost. Such a segregation has 
been considered at a scale comparable to the interpar- 
ticle distance to explain the origin of striped states in 
cupratesE'El 

In this work we will consider the opposite case in which 
the scale of the inhomogeneities is much larger than the 
interparticle distance. We consider in particular two dif- 
ferent kinds of inhomogeneities: spherical drops of one 
phase into the other phase and alternating layers of each 



phase. The first case has been pioneered by Nagaev and 
collaborators in the context of doped magnetic semicon- 
ductors in general and of manganites in particular .Bel Re- 
lated ideas have been recently presented in Refs. 

We believe that for the general understanding of the 
large scale inhomogeneous state the specific mechanism 
producing phase separation (PS) in the absence of LRC 
forces is rather unessential. Of course specific short range 
interactions in each physical system will lead to different 
A and B phases (which will also depend on the doping 
regions one considers) giving rise to different physical sit- 
uations. However, in the same spirit of the Maxwell con- 
struction, one can perform a general analysis of the phe- 
nomena due to the tendency towards PS in the presence 
of Coulomb interaction irrespectively of the microscopic 
mechanisms of PS itself. 

We consider two charged phases A and B with a com- 
pensating rigid background and we study the formation 
of inhomogeneous states in a density-driven first-order 
phase transition between A and B. By definition A and 
B have different densities; one of the phases is undercom- 
pensated and the other is overcompensated by the back- 
ground. It follows that the inhomogeneities are charged 
and they repel each other. Since the inhomogeneities are 
formed by many partides, quantum effects are negligi- 
ble and they crystallizeo The drops arrange in a Wigner 
crystal whereas the layers form a periodic structure. We 
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restrict to three dimensional (3D) textures. A large num- 
ber of small inhomogeneities minimize the Coulomb en- 
ergy but they cost too much surface energy. The distance 
between the inhomogeneities and their size are found by 
minimizing 9-ifree energy which takes into account both 
these effects. B'tl 

In ordinary PS the Maxwell construction (MC) is in- 
voked to find the range of density < n < in which 
a system prepared with the overall density n separates in 
two regions with densities n\ and vPg respectively. We 
generalize here the MC and derive the equations that 
should be satisfied in the mesoscopically inhomogeneous 
coexistence region by the chemical potential and the pres- 
sure of each phase (Sec. ||). To this end we use an ap- 
proximation in which the density within each phase is 
assumed to be constant which we name uniform density 
approximation (UDA). We solve the equations for the 
drop geometry in the simple (but general enough case) 
in which the free energy of both the A and B phases can 
be approximated by a parabola (Sec. III). 

We define a coupling constant A given by the ratio be- 
tween the the energy cost due to surface energy plus the 
LRC interaction and the energy gain in MC PS. Only 
below a critical value Ac PS is possible. Above Ac the 
system is uniform A [B) phase below (above) a critical 
density with a lattice instability close to the critical den- 
sity. 

The characteristic size of the inhomogeinites is shown 
to be of order VAZ^ with Is an electric field screening 
length. Since A is bounded by Ac (of order 1) it follows 
that the inhomogeneities are of the order of or less than 

Is- 

For small volume fraction the drop geometry is more 
stable than the layers as expected from general surface 
energy arguments. On the other hand the layered geom- 
etry being simpler serves as a ground test for approxima- 
tions. In order to validate the UDA we solve the layered 
geometry in the UDA and in the more general case in 
which the density profile can spatially vary within each 
phase. In this last case the density profile is allowed to 
adjust minimizing a free energy which is an approximate 
functional of the local density (LDA) . Both the UDA and 
the LDA are shown to give very similar results for aver- 
aged quantities (Sec. iVA). 



To illustrate the generality of theses ideas we consider 
some applications in paper II of this series. 



II. FREE ENERGY AND COEXISTENCE 
EQUATIONS: THE UNIFORM DENSITY 
APPROXIMATION 

We consider a density-driven first order phase transi- 
tion in the presence of the LRC interaction and surface 
energy. We look for the formation of a mixed state by 
increasing the density from the uniform A phase. We use 
two different geometries for the mixed state, i) The drop 
geometry consist of a Wigner crystal of drops of B phase 



in the host phase A. ii) The layered geometry is made of 
a periodic array of alternating layers of A and B phases. 

For both geometries the electronic density within each 
single phase region is taken as uniform (UDA) and in 
general it will result different from the compensating 
background density. This is of course an approximation 
since both densities will tend to adjust within each phase 
also to make the total electrochemical potential constant. 
The UDA will be relaxed in Sec. IV for the layered geom- 
etry by minimizing a free energy functional on a simple 
LDA. We anticipate here that both the UDA and the 
LDA give very similar results thus justifying our exten- 
sive use of the UDA here and in paper II. 

We start by computing the total free energy. In the 
same spirit of the MC we assume that the free energies of 
hypothetically homogeneous bulk phases are known and 
given by Fa and Fb - We define the mixing energy Em as 
the sum of the total surface energy and electrostatic en- 
ergy (computed below) . We work at a fixed total volume 
V and number of particles N. At a given temperature the 
total free energy is F = Fb{Vb,Nb)+Fa{Va, NA)+Em- 
We have to minimize this respect to Vb and Nb sub- 
ject to the conditions Vb + Va = V, Nb + Na = N. 
The volume fraction of B phase is x = Vb/V. We can 
work with the free energies per unit volume / = F/V, 
Cm = Ejn/V, Jb = Fb{Vb,Nb)/Vb and work with the 
densities ub = Nb/Vb etc. so the function to minimize 
is: 



/ = (1 - x)fA{nA) + xfBiriB) + Cr. 



(1) 



The constraint in the number of particles is written as 
n — xub + {I — x)nA and the constraint in the volume 
is satisfied by putting Vb/V — 1 — a;. It is convenient 
to define S = ub — riA and to use the constraint in the 
number of particles to eliminate hb and ha in favor of n 
and S. 

In order to compute the mixing energy we first con- 
sider the drop geometry. We assume that the drops are 
spheres of radius Rd- This will be a good approximation 
as long as x is small and the crystal field is also ap- 
proximately, spherical. This is true for fee, bcc and hep 
lattices.Ejild To compute the elee|l.rQst.atic energy we use 
the Wigner-Seitz approximation.ElE30 We divide the sys- 
tem in slightly overlapping spherical cells each one with 
the volume AttR^/S — V/Nd where Nd is the number of 
drops and i?c is the radius of the cell. Fig. |l| shows a 
schematic view of the cell density profile. 

Next we compute the electrostatic energy. The cells 
are globally neutral (by construction) and only the charge 
inside the cell contributes to the electric field in the cell. 

The charge density of phase B is ub (actually —ens 
but we drop the charge of the particles — e for simplic- 
ity). The dashed background charge density in Fig. ^ 
{—nA) compensates the A charge density ua, and a slice 
of height UA of the B charge density. For the purpose 
of computing the electrostatic energy these charge den- 
sities can be eliminated and one is left with the density 
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Ub — nA{= S) inside the drop and —{n — ua) for the 
background. We will call the former "drop contribution" 
and the latter the "background contribution". There is 
no "host" contribution due to the above cancellation. 
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FIG. 1. Schematic view of a cell density profile in the UDA 
with a drop (layer) of B phase in the host A. The origin is 
at the center of the cell. The full cell diameter (width) is 
2Rc for drops (layers). The dashed region of the background 
compensates the A density and part of the B density. 

Another assumption is that the charge is spread uni- 
formly and that microscopic discreetness effects can be 
neglected. One can see that corrections to the electro- 
static energies due to discreetness are of order a? / 
(Appendix. ^) where a is the characteristic length of the 
microscopic structure (for example a lattice constant). 
Therefore they are negligible in our analysis which con- 
siders a- 

With the above approximations the total electric field 
inside the cell is written as E = E^ + Ec; where h [d) refers 
to the background (drop) contribution. Integrating the 
square of the electric field we obtain three contributions 
to the electrostatic energy: e = + Cfc + f-d-b with 



(2) 



with eo the static dielectric constant and a similar equa- 
tion for the background. The interaction energy is 
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(3) 



The use of the static dielectric constant is well justified 
because we are assuming a static super structure which 
certainly will produce relaxation of the ions which in turn 
will screen the electric field. We are assuming by symme- 
try that the electric displacement is parallel to the elec- 
tric field. The fields can be easily evaluated with Gauss 
theorem. One obtains 
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he^Rd 
, 3 



5eoRc 
3Q^ f _J_ 
eo I 2Rc 



10i?3 



(4) 
(5) 
(6) 



where Q = —eSvd is the effective charge inside the drop. 
The volume of a drop is Vd — 47ri?j^/3 and the number 
of drops is given by Nd = VB/vd = xV/vd- We also have 
that X ~ R^^/R^. Finally the total electrostatic energy 
per unit volume can be put as: 
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i?^a;5/3(2-3a; 
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-X). 



(7) 



Setting one of the densities in 5 to zero one recovers the 
expressions obtained by Nagaev and collaborators for the 
particular case of a mixed state composed of an antifcrro- 
magnelic insulating phase and a ferromagnetic metallic 
phaseBB 

The surface energy is parametrized by a quantity a 
with dimensions of energy per unit surface. In general 
a will be a function of the densities ua, nB- The total 
surface energy per unit volume is: 



V Rc 



(8) 



These two contributions add to the mixing energy per 
unit volume e™ = Ce + Ca- 

Due to the constrain we have three parameters to de- 
termine ((5, X, Rc)- The mixing energy is the only contri- 
bution which depends explicitly on the geometry. We can 
therefore eliminate R^ in favor of 5 and x by minimizing 
Cm respect to the cell radius to get: 
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47rx(2 - 3a;i/3 + x)e^5'^ 
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(9) 



Now we consider the layered geometry. The cell con- 
sist of a layer of width 2Rc. The center of the cell is 
occupied by a layer of width 2Rd of B phase and the rest 
is occupied by A phase. Fig. |l| serves as a schematic plot 
of the density profile also in this case, r is a coordinate 
perpendicular to the layers with the origin at the cen- 
ter of the B layer. The volume fraction now is given by 
X = Rd/ Rc- By following analogous arguments as for the 
drops we obtain: 
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(10) 
(11) 
(12) 



Once Rc has been eliminated for both geometries the 
mixing energy can be put as: 
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where all the geometric information is stored in u{x) 

1/3 



w(x)=35/3(^) x(2- 3x1/3 +.x) (drops) (14) 
u{x) - (I) [3x(l - x)]2/3 (layers) (15) 

In Fig. H we plot u{x). 




FIG. 2. The function M(a;) that parameterizes the mixing 
energy for the layer geometry and the drop geometry. 

The free energy should remain invariant respect to an 
exchange of the kind A ^ B and x 1 — x. We 
will term this as "phase exchange symmetry". Fig. ^ 
shows that this symmetry is only approximately realized 
by the u{x) for drops. The deviation is due to the fact 
that the surface energy is minimized when the minor- 
ity phase inhomogeneities are spherical. Our drop so- 
lution imposes this at small x but violates this in the 
opposite case of a; ^ 1 where the minority phase inho- 
mogeneities have the complicated geometry between the 
spherical drops. In practice, however, our approxima- 
tions can be meaningful even at intermediate and large 
X because the present u{x) is approximately symmetric 
around x — 1/2. This is due to the fact that the electro- 
static energy Eq. (0) correctly cancels in this limit driving 
the total mixing energy to zero. 

A better treatment should allow at a; > 1/2 for a switch 
from the unoptimized interstitial geometry to a spheri- 
cal geometry with an energy gain given by the reflected 
u{x) at small x in u(l — x) as shown in Fig. ^ A com- 
parison between the reflected curve and the original u{x) 
shows that this geometry optimization compared with an 
apparently very bad geometry gives rise to a modest low- 
ering of the energy. The same happens when we switch 
from the layer geometry to the spherical drop geometry 



as shown in Fig. |[ We can conclude that the dependence 
on geometry is week. 

The spheric drop geometry has lower energy than the 
layered geometry as expected from general arguments on 
surface tension. The exception is close to x = 1/2 where 
our spherical drop solution is not adequate in any case. 
In fact in this region drops and the crystal potential will 
be far from spherical. The problem of establish the op- 
timum geometry close to x = 1/2 is beyond the scope of 
this work however we expect minor corrections to ther- 
modynamic quantities due to the small sensibility of u{x) 
to dramatic changes in the geometry as illustrated above. 

Although the layer solution has higher energy due to 
its simplicity it is an excellent test ground for checking 
the approximations. We take advant age o f this fact to 
test the UDA approximation in Sec. IV A. In addition 



the layer geometry has the extra advantage that, by con- 
struction, respects the phase exchange symmetry. 

Anyway since u{x) depends weakly on geometry our 
results for macroscopic thermodynamic quantities will be 
largely independent of the geometry itself. When possible 
we present our results in a geometry independent way by 
leaving the function u{x) unspecified in our expressions. 

Minimizing the free energy respect to S and x one ob- 
tains: 
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(16) 
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3(eofT)i/3 + d^J ■ 
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Here pA = -Ja + P-AnA, {p-A = dfA/duA), etc. are the 
"intrinsic" pressures (chemical potentials) of each phase. 
The word "intrinsic" stands for the values of these quan- 
tities in the presence of a fictitious fully compensating 
background, in other words they refer to a uniform single- 
phase situation. Equations (|l6|),([l^) determine the jump 
in these quantities at the interface in order to have ther- 
modynamic equilibrium when long range Coulomb forces 
and surface energy are present. These equations are valid 
for a general geometry described by the function u{x). 
Notice that as long as u{x) preserves the phase exchange 
symmetry the equations also preserve this symmetry. 

To analyze the effect of the long range forces and 
of the surface energy in the jumps let us neglect for 
simplicity the density dependence of the surface energy 
{da/duA — da/duB — 0) and concentrate on the drop 
geometry. Due to the different charge distributions, 
the electrostatic potential energy —ecf) of an electron in- 
side and outside the drops is different. In equilibrium 
this jump in the electrostatic potential should be com- 
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pensated by a jump of the intrinsic chemical potentials 
[Eq. to make the electrochemical potential constant 
in the whole system. For S > the drop repels electrons 
so the electrostatic potential energy will be lower outside 
the drop i.e. —ecfiA < —eft's- The intrinsic chemical po- 
tential outside will have to be larger than inside as the 
sign in Eq. ( [l6| ) implies. 

Regarding the pressure, in equilibrium the intrinsic 
pressure inside the drop, pB, should equal the pressure 
exerted by the host pA plus the pressure exerted by the 
mixing forces. For (5 > the electrostatic energy induces 
a negative contribution to the pressure since an increase 
in the drop volume at constant particle number decreases 
the difference in densities between the interior and the 
exterior of the drop and hence the Coulomb cost. This 
effect is given by the first term in Eq. (pT|). The second 
term proportional to u'(x) is a geometric contribution. 
Both terms are discussed in more detail in a specific ex- 
ample in Appendix ^ 

In the limit e ^ one gets ^.b — l^A ~ M and pA = 
Pb — P i-C. pl5 = Jb — J A which arc the conditions for 
MC. 



III. GENERAL ANALYSIS OF THE MIXED 
STATE IN THE UNIFORM DENSITY 
APPROXIMATION 

In this section we set up the basic ideas for inhomoge- 
neous solutions. For simplicity we model each phase free 
energy with a parabola and we assume that the surface 
tension is density independent. Without loss of gener- 
ality we write the parabolas as a quadratic expansion 
around the MC densities: 

1 



(18) 
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The quantities with the "0" superscript (or subscript 
below) satisfy MC in the absence of LRC forces i.e. 
/s ~ — M°<^o and (5o = rig — n^- The linear sloop 
is the same for the two phases due to the MC con- 
dition. The MC density no and the volume fraction are 
related by uq — -I- Sqx. The constants kA, k^ are 
essentially the compressibilities of the two phases.Ej 

For non interacting electrons at T=0 the compressibil- 
ity coincide with the density of states. For the 3D free 
electron gas we have: 
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'^free 



(19) 



with m the electronic mass. 

Another useful realization is a nondegenerate gas 
where we have: 
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(20) 



Our aim in the following is to obtain the equations 
which control the deviation from MC behavior in the 
presence of the mixing energy. 

We define a dimensionless global density 

n = {n- n\)/So 

which measures the distance from the point in which B 
phase appears in the absence of Coulomb forces. In MC 
the coexistence region is given by < n' < 1. 

Eqs. (p^,([l7|) determine 6, and x for a fixed density 
where now fiA, (J-b^ Pa, and can be expressed in terms 
of the parameters appearing in Eqs. (|l8|). 

In practice it is much easier to solve the equations by 
fixing the volume fraction x and solving for S, and n, 
i.e. we find which density one should put in the system 
to obtain a mixed state with a given volume fraction. 
This is because the solutions happen to be multivalued 
functions of n whereas they are single valued functions 
of X (see below). 

For a fixed volume fraction x we define the dimen- 
sionless density deviations from the MC values: n = 
{n — no)/(5o and S = {S — (5o)/(5o. The density devia- 
tion h measures the shift in the global density needed to 
have the same volume fraction of a system without LRC 
interaction. 

To fix the energy units it is convenient to choose one of 
the two compressibilities as a reference, for example the 
largest. We define km — max(fcyi, /c^). Energies per unit 
volume will be measured in units of the characteristic 
PS energy S^/km- The latter is essentially the difference 
between the uniform parabolic free energy and the MC 
free energy at the characteristic density Sq. 

Now we define two important reference lengths scales 
in the theory. The characteristic size of an inhomogeneity 
for which the Coulomb energy balance the surface energy 
is given by the Rc of previous section with the geometric 
factors dropped and the density evaluated at the MC 
value. This define the scale 
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(21) 



The other length is given by = eo / {Ane'^ k,n) ■ By re- 



laxing the UDA we will show in section Sec. IV that 1^ is 
a screening length. In other words if the reference phase 
(the one with km) is interpreted as a metal Is is the char- 
acteristic distance in which the electric field penetrates. 

The theory has two dimensionless parameters. One is 
the ratio /cb/A:^. The other measures the strength of the 
mixing energy energy effects in units of the characteristic 
PS energy 5'^ /km and is given by: 
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o(^) (?) (22) 



5 



The characteristic mixing energy is given by the factor 
with the power 1/3 in the middle expression. The con- 
stant A characterizes the competition of the mixing en- 
ergy cost and the MC Uke energy gain due to phase sepa- 
ration. The couphng constant goes to zero as e with 
a finite. This correspond to the usual PS. The case cr ^ 
with finite e correspond to an infinite number of drops 
(or layers) of zero radius. In this maximum intermixing 
situation the charges of the two phases spatially coincide 
and the Coulomb cost goes also to zero so that the MC is 
again valid. Notice however that this last idealized situa- 
tion cannot be reached in practice because at some point 
for small drop radius the continuous approximation will 
fail. 

Inserting the explicit expressions [Eqs. (^sj)] of Ja and 
Jb in Eqs. (|l6|),(|l^) we obtain the following equations 
for the density deviations: 
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Eqs. ( p3| ) can be solved numerically for general values 
of the parameters. For small A i.e. for small mixing 
energy we can linearize the equations neglecting higher 
order terms in 5 and h. We will refer to this as the 
linearized UDA. We get: 
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For the sake of simplicity in the following we will con- 
sider the linearized solution. We checked that for all the 
physical properties the difference between the linearized 
and the exact solution is quite small in the range of A 
where the drop solution is stable. 

The linearized solution takes a simple form and is ex- 
plicitly symmetric respect to an exchange of phases when 
written in the original variables: 
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In the case of A = 0, according to MC, the system sep- 
arates in two phases with densities n^^, rt^ respectively 
independently of the volume fraction. For nonzero A and 
small X the B phase divides in drops or layers and the 
density in each phase depends on the volume fraction of 
B phase. The deviation of each density from MC pre- 
diction is proportional to A and to the compressibility of 
each phase. Notice that the density of an incompressible 
phase (fc 0) does not depend on the volume fraction 
even in the presence of LRC forces. 

In Fig. ^ we show the behavior of the two functions 
which determine the dependence of the densities on the 
volume fraction. In the drop geometry and for small x 
both ua and ub tend to be larger than in the MC case 
whereas in the layered geometry only ua is larger. This 
gives rise to minor qualitative differences in the behavior 
of drops and layers. Apart from this the overall behavior 
is similar. 
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FIG. 3. The dimensionless functions that determine the 
change in ua (upper curves) and ns (lower curves) for small 
A vs. the volume fraction x for the layered and the drop 
geometry. [See Eq. ^)] 



The equation for the density of one phase [Eq. (|2£ 
has a transparent interpretation in the limit in which the 
other phase, say A, is incompressible (fc^ = 0). This case 
is solved in detail in Appendix Assume that A phase 
is the vacuum and so exerts no pressure and has zero 
density. We can consider that the mixing forces due to 
the electrostatic and surface energies exert an "external" 
pressure on the B phase inhomogeneity. In equilibrium 
the intrinsic pressure of i? phase (pb) should compensate 
this "mixing pressure" [pB — Pm)- The latter is shown 
in Appendix M to be given by: 



6 



dx 



2p 
3a: 



(27) 



On the other hand a change in the external pressure cor- 
respond to a change in the density according to the 
B phase equation of state. TJiiis follows directly from our 
definition of compressibility 113 

, _ ^"-B 
«B = rig- 

where we have replaced a derivative by a finite different 
ratio. We can obtain the second linearized expression in 
Eq. (^ ) directly from this definition using that the MC 
density correspond to zero intrinsic pressure. 
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The mixing pressure can be negative as explained in 
the Appendix This implies that the density is less 
that the MC density. From the lower cures in Fig. ^ we 
see that for drops the pressure is positive for small x and 
then becomes negative whereas for layers the pressure is 
negative for all x. 

Remarkably in both cases the mixing pressure is a de- 
creasing function of x. Since in general x is an increasing 
function of n' we can anticipate that n b will decrease as 
n' increase (see below). Notice that for small x we have 
Pb ~ ii'(a;)/3 so a decreasing mixing pressure can be di- 
rectly related to the negative curvature of u{x) (Fig. ||). 

Coming back to the general solution in Eq. (|2^) we 
are interested in the dependence of these quantities as a 
function of the global density n' , our true control vari- 
able, rather than as a function of the volume fraction. 
Hence we need the volume fraction as a function of the 
global density n' . From the solution of the linearized 
equations we find: 
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(28) 
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Since all physical quantities depend on the densities this 
completes the solution of the problem. 

Specific results will be presente d in t he next section 
for the drop geometry and in Sec. IV A for the layered 
geometry. 



A. Results of the UDA for the drop geometry 

Now we consider the drop geometry and we analyze in 
detail the two cases: i) the compressibilities of the two 
phases are equal (fc^ — ~ fcm) and ii) one of the 
compressibilities is zero. 

In Fig. 1^ we plot the volume fraction as a function of 
global density from Eq. ( p8| ) for the drop solution. The 



volume fraction is a multivalued function of n' and in 
the case ks = kA has a lower branch close to a; = 0, 
an intermediate branch, and an upper branch close to 
X — 1. The intermediate branch is the physical solution. 
This will be shown below by looking at the free energy. 
The physical solution has the intuitive property that the 
volume fraction increases as global density increases. 
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FIG. 4. Top panel: Volume fraction vs. n' for (from left to 
right at the bottom) A = 0, 0.1, 0.2, 0.3, 0.4, 0.5 and ks = kA- 
For A = 0.4 we indicate with a vertical line the discontinuity 
in the volume fraction to go from the uniform solution to the 
drop solution by increasing the density. Bottom panel: Same 
for ks =0. The approximations done are rigorously valid 
only for small x. 

We see that the bifurcation density nj^j^ at which the 
phase separated solution appears for finite A is larger 
than in MC. On the other hand B phase appears with a 
finite volume fraction and its growing rate is larger than 
in the MC case. Remarkably both the volume fraction 
at the bifurcation point and the bifurcation density 
are almost the same for ks = kA and for /cs = 0. They 
depend only on A as can be seen by comparing the two 
panels in Fig. IJ. 

In the case ks = the constraint between the volume 
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fraction and the densities together with the fact that the 
B density is fixed make all the curves to converge to the 
MC case when a; ^ 1 as shown in Fig. |[ The same 
happens when fc^i — and a; — > 0. 
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FIG. 5. fA-f, fB-f, and/-/" in the drop solution for 
A = 0.1,0.2,0.3,0.4,0.5 (from bottom to top) and ks = kA 
vs. n'. Here /° is the MC free energy for A = (a straight 
line). The cross indicates the value with a; = of the drop 
solution for each A. The black dot indicates the bifurcation 
point in which the drop solution first appears when density 
increases. 



given by ribi / vs. A as shown in Fig. ^. The uniform-drop 
boundary line is determined by the condition dn! jdx — 
(see Fig. I). 

For A > Ac the homogeneous solution is stable for any 
global density. The uniform A-B boundary line is deter- 
mined by the crossings of the parabolas in Fig. ||. 

When one of the compressibilities goes to zero, say fc^, 
the crossing moves to the right in Fig. |^ and the uniform 
B region shrinks until the boundary line for uniform B 
phase approaches the MC value (n = n^). At the same 
time Ac increases. Analogous results are obtain for kA 
going to zero. 
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To decide the stability of the solution we have to 
compare the drop solution with the single phase solu- 
tion. In Fig. 1^ we show /^(n'), jsin') and the total 
free energy with ks = kA for various A. The MC line 
/O(n') = fA+n'ife - Ia) ^las been subtracted. The en- 
ergy also is a multivalued function of n' . As the density 
increases the drop solution appears at n'^jj (indicated in 
Fig. H by a black dot) with two different branches. In 
the upper (unstable) branch x decreases with density till 
the point x — Q highlighted with a cross in Fig. |[ For 
the lower branch one finds the expected behavior i.e. x 
increases with density. The upper branch is almost de- 
generate with the bulk fA{n-) free energy. Near the bifur- 
cation the three solutions (homogeneous, drop stable and 
drop unstable) are very close in energy. Approximation in 
the solution of the Eq. (|3|) can lead to wrong conclusions 
about the relative stability. In this case one has to refer 
to the non-linearized solution. For the latter (not shown) 
we find that the bifurcation density ni,i / is lower than the 
density Uc at which the energy of the lower energy drop 
solution crosses the energy of the uniform phase /yi(n'). 
However the difference between ric and Ubif is negligi- 
ble for all practical proposes except for the largest A. In 
this case there is a small region (Ac = 0.49 < A < 0.57) 
in which the lower energy drop solution still exist but is 
less stable than the homogeneous solution. If we neglect 
this small effect the phase diagram of the drop solution is 



FIG. 6. Locus of existence of the low-energy drop solution 
in the \-n' plane for kA = ks- This almost coincides with the 
phase digram in the sense that when the drop solution exist 
it is more stable than the uniform solution except close to 
A = 0.5 and in a very narrow region around the drop-uniform 
boundary line (see text). 

In the upper panel of Fig. we show the density of each 
phase as a function of the global densities for ks = kA- 
Increasing the global density the transition occurs from 
the uniform A phase, with density higher that the MC 
one, to the drop state. In the MC case the density of 
A phase is continuous at the transition and remains con- 
stant in the coexistence region. For nonzero A the A den- 
sity has a discontinuity when the drops occur (see inset). 
Remarkably both local densities decrease as the global 
density increases. This is due to the behavior of the mix- 
ing pressure as explained above and in Appendix |b[ In 
the case kA — the regions whit ub > (ns < n^) 
can be directly associated with positive (negative) mixing 
pressures. 

Compared to the upper panel the lower curves for nA 
shrink to the MC case and the upper curves for hb remain 
very similar (even quantitatively) except close to n' 0. 
We mention that in the case fc^ = (not shown) a similar 
effect is seen exchanging A with B. 

In Fig. H we show the cell radius and drop radius in 



8 



units of the screening length as a function of density for 
ks = kA- Both the cell and the drop radius are typically 
on the scale of a few screening lengths Is for not too small 
A and have a finite size at the appearance of the mixed 
state. The cell radius decreases as the density goes away 
from the bifurcation value to reach a minimum close to 
n' = 1/2. The minimum would be exactly at n' — 1/2 in 
an exact computation due to phase exchange symmetry. 
This is show below in the layered solution. 

The (B phase) drop radius instead is intrinsically 
asymmetric and increase monotonously with the density 
reflecting the transformation of the cell from A phase to 
B phase. 

For A ^ the cell radius and the drop radius behave 
as i? ~ VXls ^ [creo/((5oe)^]^/^. As stated in Sec. ||they 
diverge as e ^ indicating that MC can be realized with 
a single large drop of B phase in A. 
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FIG. 7. Normalized densities of each phase as a function 
of normalized global density n' for different A. The upper 
panel is for = kA and the lower panel is for fc^ = 0. For 
each panel the lower curves correspond to A phase and the 
upper curves to B phase. In the coexistence region multi- 
valued densities appear. The long branch is the physical one 
and the short branches are unphysical. The inset shows an 
enlargement of the A density to resolve the discontinuities. 



Another peculiarity of the curves in Fig. g is that the 
free energy of the drop solution has the "wrong" curva- 
ture, that is the compressibility (defined from d'^f/d'^n) 
is negative. This does not necessarily imply an instability 
since the usual stability condition of positive compress- 
ibility is formulated for a neutral system, that is including 
the background compressibility. Since we are assuming 
the inverse background compressibility to be an infinite 
positive number (in our analysis the background density 
has a fixed homogeneous value) it follows that the total 
compressibility is positive and from this point of view the 
system is in a stable mixed state. Of course this does not 
guarantee stability against more complicated solutions 
than the simple crystal of drops. 




0.2 0.4 0.6 0.1 
n' 




FIG. 8. Rc and Rd in units of the screening length Is de- 
fined above Eq. (^2|) vs. n' for fcs = kA- We show the curves 
for A = 0.1,0.2,0.3,0.4,0.5 which increases from bottom to 
top in the top panel and from right to left at the top in the 
lower panel. In the top (bottom) panel for each curve the 
lower (upper) branch is the stable one. 

The situation is more severe for A > Ac where the drop 
solution, if it exists, is not stable. In this case, the sys- 
tem remains always single phase and the free energy is 
given by the branches of the parabola with the smaller 
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energy in Fig. ^. It changes suddenly from A phase to 
B phase at the density n'^ for which fA{n'^) = fB{n'c)- 
(For our parameters = 0.5). The problem is now that 
the energy has a cusp pointing upwards at n'^ which im- 
plies an infinite negative inverse compressibility. This 
will compete with the infinite positive inverse compress- 
ibility of the background. Clearly one should consider 
in this case the background compressibility (e.g. the lat- 
tice compressibility) since the beginning. As a first step 
we can add to the above electronic free energy a back- 
ground free energy contribution fb{n) = (n — /2ki,. 
A very rigid (but not infinitely rigid) background is de- 
scribed by a very small fcb > which correspond to a very 
narrow parabola for the background free energy. The to- 
tal free energy, background plus cusp, will have a cusp 
pointing up with two local minima nearby. Since now 
the total free energy corresponds to a truly neutral sys- 
tem one can make a MC between the two local minima. 
One obtains a phase separation between A and B with 
the background adjusting its density in each region to 
the density of each phase to make it neutral. The same 
argument applies at the critical density where the drop 
solution crosses the uniform solution, although the neg- 
ative dip is much less pronounced in that case. Usually 
the electronic system is a crystal where the background 
is provided by the ionic lattice. If one trays to prepare 
the crystal with an electronic density close to the critical 
one the system can break in two pieces each one with a 
different lattice constant. Typically the crystal is not at 
a fixed volume but at a fixed external pressure P. (We 
use capital P to distinguish the pressure exerted on the 
crystal as a whole from the electronic pressures of the 
phases pA and pb)- In this situation MC determines the 
equilibrium pressure Pq for phase coexistence. Pq will 
depend on the global doping so above Ac, Pq vs. doping 
determines a phase boundary line which will cut ambient 
pressure at some critical doping. 

Since the electronic free energies depend on external 
parameters, a remarkable implication is that the critical 
doping will also depend on external parameters like mag- 
netic field, temperature, pressure, etc. In other words a 
crystal can be driven from a single phase to a two phase 
situation by changing external parameters. This is very 
reminiscent of the situation in some manganites where 
one finds that a single phase crystal brakes in a multido- 
main crystal by lowering the temperature. The multido- 
main system shflssjs lattice mismatch and large stress at 
the interfaces. E3o 



In Sec. [VA analogous results are presented for the 
layered geometry case and compared with a more elabo- 
rate computation which relax the UDA. In paper II we 
apply to different physical systems the ideas developed 
in this section. 



IV. LOCAL DENSITY APPROXIMATION 

In this section we generalize our results to take into 
account the full spatial dependence of the density. The 
basic assumption is that we can write the free energy of 
each phase as the spatial integral of a free energy density 
which is a function of the local density, i.e. we are using 
a local density approximation (LDA). The free energy 
reads: 

F = [ d?V,fAHv)]+ [ dPvfBHv)] 
JreA JreB 
1 



d^rE^ -f- Sabct 



(29) 



Here r e A indicates that the integral is restricted to the 
regions of phase A and Sab is the total interface surface 
between A and B and we assume for simplicity eg = 1. 
One should be careful not to double count in a surface 
energy costs that are due to the spatial variation of the 
charge since this will be explicitly taken into account in 
the first three terms. On the other hand one can include 
in a other effects, like magnetic, which would not be in- 
cluded otherwise. For simplicity we will assume a to be 
density independent. 

The electric field is related to the total charge density 
(electronic plus background) through Poisson equation: 



V.E = Anp 
with the total charge density: 

p ~ — e[n(r) — 



(30) 



(31) 



Here n is the global density of the previous section and 
the bar distinguishes it from the spatially varying den- 
sity n(r). Notice that en is the charge density of the 
background. The condition of neutrality is written as: 



1 
V 



reA 



d^rn(r) 



1 

V 



d^rn{r) 



(32) 



Using 71 (r) — ua for r £ A and ?i(r) = tib for r e P one 
recovers the UDA. 

Instead of minimizing the functional with respect to 
the density it is convenient to use Eqs. (|30|),(^ to 
express the density as a function of the electric field, 
[n = n(V.E)] and minimize the functional with respect 
to the electric field profile. We look for periodic solutions 
(layer, crystal, etc) and restrict the computation to one 
cell. 

Minimizing the free energy [Eq. (|2|)] respect to the 
electric field one obtains: 



E = -iv^[n(V.E) 
e on 



(33) 



where X — A or B when r £ A oi r G B respectively. 
This differential equation together with the boundary 
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condition determines the field profile. The boundary con- 
dition at the cell boundary and at the internal boundary 
will be discuss in the example below. Once the electric 
field profile is known for a given geometry the density 
profile is given by Poisson equation. As a final step one 
should optimize the geometry. 

Introducing the parabolic expressions [Eqs. (Ilq )] to pa- 
rameterize the free energy densities in Eq. (P3[) one ob- 
tains: 



E = ZjfVV.E 



(34) 



with — {Ane^kx) ^- Clearly is the screening length 
as anticipated in Sec. [II. If we use the compressibility 



of a free electron gas for kx [Eq. |19| and reintroduce the 
dielectric constant Ix corresponds to the Thomas-Fermi 
screening length: 



/2 - 



1/3 



4e2m(n'^)i/3 ' 



(35) 



We reach Thomas-Fermi theory which is the simplest ver- 
sion of the LDA used for electronic structure computa- 
tions. If we use the nondegenerate gas compressibility 
[Eq. (EOh] Ix is the Debye-Hiickel screening length. 



(37) 



PB 



Eq cosh{r/lB) 
AttIb smh{Ra/lB) 



The electric field at the A-B boundary can be related 
to the jump in the density at the interface: 



-47re[nB(i?d) - riAiRd)] 



[Ib tanh(i?rf/;B)]-i + {Ia tanh[(i?c 



Rd)/lA]}-' 
(38) 



It plays the same role as the parameter S in Sec. || so that 
we can find the optimum charge distribution between A 
and B by minimizing the free energy with respect to Eq. 

After replacing Eqs. (37),(^) in the expression for the 
free energy [Eq. (p9[)] and minimizing respect to Eq we 
find: 



Eq = 



AttcSq [Ib^ [n' - 1) - lA^n'] 
Ib/ tanh(xi?c/^s) + Ia/ tanh[(l — x)/li 



(39) 



where 5q and n' are defined as in Sec. |l| and Rd has been 
eliminated in favor of the volume fraction with i?^ = xR^. 



A. Solution for the layered geometry 

In the layered geometry the differential Eq. ( ^^ re- 
duces to a one-dimensional problem and can be readily 
solved. The geometry is identical as in the UDA approx- 
imation (Fig. |l|). The central B layer has width 2R^ and 
the cell has width 2i?c. The r coordinate is perpendic- 
ular to the layers and r — Q corresponds to the center 
of the B layer. By symmetry the field is zero at r = 
and r — Rc. In this case the boundary condition — 
for the electric field perpendicular to the surface at the 
cell boundary automatically warrants the neutrality con- 
dition [Eq. (^2|)] due to Gauss theorem. 

Apart from the cell boundary the cell itself has an in- 
ternal boundary that divides A and B phases we call Eq 
the electric field at the A-B boundary. The value of Eq is 
also optimized and this provides an additional boundary 
condition. 

The solution is of the form: 



EA{r) 
EB{r) 



Eq 



Eq 



sinh[(r - Rc)/Ia] 
sinh[(i?<i - Ra)/lA 

sm]i{r /Ib) 
sm\\{Rd/lB) 



(36) 



where EA{r) = E{r) for r € A, etc. 
The charge density is given by: 



PA 



Eq cosh[{r~R,)/lA] 
AttIa smh[{Rd - Rc)/Ia] 




FIG. 9. Volume fraction vs. n' for (from left to right at 
the bottom) A = 0,0.1,0.2,0.3,0.4,0.5 and fcs = /ca in the 
LDA (thick line) and the UDA (thin line). Only the lower 
left corner of the plot is shown since the upper right corner is 
symmetric by phase exchange. For the U DA ap proximation 
the lower branch is unphysical like in Sec. [II A 



At this point the total free energy per unit volume 
/ = F/V takes the form: 



(T 



(40) 



+ 2TTe^dlfA{n'f{l -x) + llx{l ~ v!f\ 

i?c{ZB/tanh(a;i?c/ZB) + ;A/tanh[(l - x)Rc/lj^} 
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The first two terms are the MC free energy, the third 
term is the surface energy and the last two terms are both 
contributions due to the shift from the MC densities and 
due to the electrostatic energy. 

The last step is to minimize this free energy with re- 
spect to the volume fraction and Rc. This gives two 
equations which can be solved numerically for Rc and x. 
As in Sec. Ill it is easier to fix x and solve for Rc and n' . 

In the following we present results for the case ks — kA 



and compare with the linearized UDA of Sec. Ill for the 
layered geometry. 




FIG. 10. /a ~ f, fs ~ f, and / - /" in the layered 
solution for A = 0.1, 0.2, 0.3, 0.4, 0.5 (from bottom to top) and 
ks = kA vs. n' in the LDA (thick line) and the UDA (thin 
line). Here is the MC free energy for A = (a straight 
line) . 

In Fig. H we plot the volume fraction as a function of 
global density in the LDA approximation and the UDA 
approximation. Clearly the result are very similar even 
quantitatively. In the UDA there is a jump on the vol- 
ume fraction form zero to a finite value. In the LDA 
the volume fraction is not discontinuous but grows very 
rapidly at the threshold for the appearance of the inho- 
mogeneous state. Another important difference is that 
the solutions are not any more multivalued in the LDA. 
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FIG. 11. Normalized spatially averaged densities of each 
phase as a function of normalized global density n' for difler- 
ent A, fcs = kA and the linearized UDA (thin lines) and the 
LDA (thick lines). The lower curves correspond to A phase 
and the upper curves to B phase. In the coexistence region 
multivalued densities appear in the linearized UDA. The long 
branch is the physical one. 

In Fig. |l^ we show fA{n'), /s(n') and the total free 
energy with fc^ = kA for various A. The MC line 
fin') = fA+ "'(/s - /a) has been subtracted. The 
behavior of the layered solution in the U DA is similar 
to the one found for drops in Sec. Ill A and coincides 
with it at small A. In the LDA approximation multival- 
uation disappears. The relaxation of the UDA approx- 
imation produces obviously a gain in energy since the 
functional that we are minimizing is the same in LDA 
and the UDA but in the UDA we are imposing an extra 
constrain on the densities. The gain in energy however 
is quite small. The phase diagram in the UDA and the 
LDA (not shown) are both very s imilar (even quantita- 
tively) to the one for drops of Sec. Ill A except that they 
are fully symmetric. The critical A above which the inho- 
mogcneous solution is never stable is given for ks = kA 
by Ac = (9/5)1/3/2 - 0.61 in the LDA and by Ac - 0.70 
in the UDA. 

In Fig. ^ we show the densities in each phase in the 
UDA. This is compared with the densities of each phase 
in the LDA averaged spatially over the space spanned 
by each phase. Again the behavior is remarkably similar 
and the density discontinuities of the UDA become very 
steep changes with LDA. 

Finally in Fig. [ij we show the behavior of the dimen- 
sions of the cell and of the B layer as a function of global 
density. Due to perfect phase exchange symmetry the 
cell width 2Rc as a function of n' is symmetric and has 
the minimum exactly at n' = 0.5. The discontinuous 
jump at the threshold in the UDA becomes a divergence 
in the LDA. For the same parameters the cell width are 
smaller in the UDA than in the LDA. This can be un- 
derstood by noticing that in the UDA the widths are of 



12 



order Ig = [a/{5oe)'^Y^^ . Roughly speaking we can say 
that the effect of the LDA is: i) to increase the surface 
energy due to the bending of the charge distributions at 
the surface and ii) to screen the electric fields which can 
be schematized as an effective reduction of the charge e. 
Both effects tend to increase the with of the layers as 
found. 

For small A Fig. |l| shows that the LDA and UDA ra- 
dius coincide just as the full solution. This is because 
Id ~ VXls « Is [c.f. Eqs. so that the den- 

sity is almost constant inside the layer even in the LDA 
and the solutions are virtually the same. In this case the 
Thomas-Fermi approximation is ineffective to generate a 
surface energy since all surface energy effects other that 
the ones explicitly included in a are due to density varia- 
tions. In other words if one sets a — the system prefers 
to make small drops to avoid both the Thomas- Fermi in- 
duced surface energy effect and the Coulomb cost. This 
however is a drawback of the Thomas- Fermi approxima- 
tion since small drops will certainly have a large surface 
energy due to the confinement of the electron gas. It is 
well known that Thomas- librmi theory is a poor approx- 
imation to model surfacesJl3 

If one increases A inhomogeneities are possible until 
the point in which Id ^ Is and A = Ac. It is not possible 
to have inhomogeneities of dimension Id » Is because 
in the region far from the surface screening makes the 
local density to coincide with the global density and this 
inhibits any PS energy gain. It is then convenient for the 
system to avoid any surface and remain single phase. 
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FIG. 12. Rc and Ra in units of the screening length Is 
defined above Eq. ( |2^ ) vs. n'. in the linearized UDA (thin 
Unes) and the LDA (thick lines). We show the curves for 
A = 0.1,0.3,0.5 which increases from bottom to top in the 
top panel and from right to left at the top in the lower panel. 

In Fig. |l^ we show the density profile for A = 0.3 
and for two different values of the global density. One 
is close to the threshold for the appearance of B phase 
(n' = 0.353). In this case the A density is close to the 
density of the background and bends down close to the 
interface to screen the B layer charge. Well in the bulk 
of A phase, where the charge density coincides with the 
density of the background, we have E ^ as expected for 
a metal. When the global density increases the local den- 
sities decrease according to the behavior discussed before 
for the average densities (Fig. |l^). The layers become of 
the order of the screening length and the electric field is 
never completely screened. 
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FIG. 13. Density profile for A = 0.3 and different values 
of the global density. The region close to r = corresponds 
to B phase and the rest is A phase. The structure repeats 
periodically in the r direction. The horizontal lines signal the 
global density. 



V. CONCLUSIONS 

In this work we have generalized the Maxwell construc- 
tion to a situation that appears often in strongly corre- 
lated electronic systems, i.e. phase separation frustrated 
by the LRC interaction. 

We discussed i) the stabilization of the uniform phases 
as the frustrating forces are increased, ii) the anomalous 
behavior of the frustrated phase separated mesoscopic 
state and iii) the singular behavior which results in a 
lattice instability when frustration dominates. 

We used a UDA and a more involved LDA approach. 
Both are shown to give very similar results thus justifying 
in general the use of the much simpler UDA. For the LDA 
we have approximated the energy functional in the case 
of a metal with the the simplest LDA functional i.e. the 
Thomas-Fermi approximation. Our formulation however 
is general and allows for more sophisticated functionals. 

As it is intuitively expected, the LRC interaction tends 
to stabilize the non-separated uniform phases in the pres- 
ence of a rigid background. This has been illustrated in 
the general analysis of two generic phases described by 
parabolic free energies. We have shown that the region 
of phase separation contracts when the LRC and surface 
energy effects are switched on and disappears above a 
critical value of a dimensionless parameter A. This pa- 
rameter plays the role of an effective coupling and char- 
acterizes the competition between the energy cost due to 
the surface and Culombic energy and the energy gain in 
the MC i.e. it controls the degree of frustration. The bal- 
ance between these energies determine whether the phase 
separated state exists or not. 



When A is small (A < Ac) a mixed state arises. We have 
modeled this situation by considering a Wigner crystal of 
drops of one phase hosted by the other phase and a lay- 
ered geometry which behaves as one dimensional analog 
of the Wigner crystal. We believe that our general con- 
clusions (including the existence of a critical A) are not 
sensitive to the geometry of the mixed state as long as 
the two length scales Rc and Rd are present and both 
are much larger than the interparticle distance. The for- 
mer length (cell size) characterizes a periodic structure 
and the latter (bubble size) how this periodic structure 
is divided to host the two phases. An indication that 
the geometry is not very important comes from the fact 
that the plots of the physical quantities in Sec. Ill for 
kA = are quite symmetric to an exchange of the two 
phases, each one having a different shape. This means 
that the behavior of the drops is not much different from 
the behavior of their counterpart, the interstitial regions. 
The same happens when comparing the behavior of drops 
and layers. 

In the mixed state novel non-linear effects appear 
which are not present in the unfrustrated MC. Within 
the UDA the volume fraction and the drop radius of the 
minority phase do not start from zero but from a finite 
value and the transition to the drops state is abrupt. In 
the LDA physical quantities are not discontinuous but 
grow very steep at the threshold mimicking the discon- 
tinuous behavior. 

A further non-linear effect in the drop state is that 
the local densities of each phase have an anomalous be- 
havior decreasing as the global density increases. This 
can affect properties of the system which are sensitive to 
the local density and will be illustrated in paper II with 
the Curie temperature of the manganites. We empha- 
size that also local probes like NMR, core spectroscopy 
etc. should be sensitive to this effect and may be used to 
detect Coulomb frustrated phase separation in real sys- 
tems. 

In the case of strong Coulomb interaction and large 
surface energy (A > Ac) a transition between two uni- 
form phases occurs. We have shown that in this case 
the compressibility is singular and a lattice instability 
will take place if the ionic background is not fully rigid. 
The system (both electrons and ions) can separate in two 
neutral phases with different specific volumes. 

In principle also at the transition point to the drop 
state a lattice instability can arise for the same reasons 
discussed in the A > Ac case, although the instability is 
now much weaker. 

When do we expect such a mesoscopic phase separa- 
tion to prevail against microscopic phase separation (like 
stripes)?. In order to have mesoscopic phase separa- 
tion we need that the interparticle distance (~ n~^^^) be 
smaller than the inhomogeneous length which should 
be smaller than the screening length Is- This implies: 
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We see that large dielectric constants favor both large 
drops and inhomogeneous states so polar materials which 
have typically large static dielectric constants (eq ^ 10 
100) are ideal candidates. Small Sq or large a favors large 
drops but a too small Sq or a too large a can inhibit 
phase separation at all. Small values of Sq can occurs in 
manganites where typically a variety of different ground 
states with close densities compete with each other (see 
paper II for a specific example) . This suggest either large 
drops or total frustration with lattice instabilities close to 
the transition from one phase to the other. We mention 
that these lattice instabilities, which also involve volume 
variations, are reminiscent of the macMscopic phase sep- 
aration observed in some manganitesJl3 

Finally small compressibilities favor the PS states. 
This suggest that these effects can be important for bad 
metals or close to metal insulator transitions. 

We believe that to some extent at least some of the ef- 
fects found here can survive also in the microscopic frus- 
trated phase separation. In fact for example, the correc- 
tions to the electrostatic energy due to the discreetness of 
the charge which are computed in the Appendix ^ for 
the very unfavorable case of a classical Wigner crystal 
can be irrelevant for small metallic inhomogeneities due 
to quantum blurring. In this case however one should 
take into account the structure of the underling atomic 
potential. Of course if quantum blurring effects are two 
strong one should be concern with the stability of the 
whole superstructure against quantum fluctuations. 



APPENDIX A: CORRECTION DUE TO 
DISCREETNESS OF THE CHARGE 

In order to compute the electrostatic energy in the 
UDA [Eq. (^] we assume that the charge within one 
drop is spread uniformly. Variations of density can a rise 
because of screening effects as discussed in Sec. [V A and 
because of intrinsic charge inhomogeneities internal to 
the particular phase. Here we discuss the latter effect. 

Let us now consider an extreme limit and assume that 
both phases A and B are two classical Wigner crystals 
of electrons as a prototypical case in which the charge is 
intrinsically non-uniform. What is the correction to the 
Eq. (0)?. 

In the host phase we neglected the interaction between 
the neutral A Wigner crystal of electrons and a back- 
ground of charge density (n — nA)e (see Fig. |l|). The 
fluctuation of the charge inside the crystalline Wigner- 
Seitz cell can make this interaction nonzero. Also for 
the phase forming the drop we have to consider the in- 
teraction between the neutral B Wigner crystal and a 
background of charge density {n — nB)e. 

The electrostatic contribution per drop is: 



with 



lOrA 



Na 



4l" 3. N 

lA = ~ "^''^ 



and a similar expressions for B phase. Here ^ 



T'A- ^/^A 

and Na is the number of electrons of A phase in a drop: 

1 



Na = riAVdi 1) 



Ne 



The total contribution per unit volume to the electro- 
static energy is: 



[(ns - n)nBr%x + {tia - n)nAr\{l - x)] 
o 

(Al) 

Clearly [c.f. Eq. (0)] the correction Aee/eg is of order 
r\ b/ so it is negligible unless the volume of the drop 
is of the order of the volume per particle in which case 
the whole computation has no sense. 



APPENDIX B: "METALLIC" DROPS IN 
"VACUUM" 

In this appendix we discuss in detail the case of a com- 
pressible phase {B) growing in an incompressible phase 
[kA — 0). This simplifies the physics because the A den- 
sity is fixed so there is not interchange of particles and 
the B density is not any more bivaluated for small den- 
sities as shown in the lower panel of Fig. |^. We present 
an alternative treatment of the frustrated phase separa- 
tion phenomena which enlightens the underlying physics 
and discuss the pressure exerted by the mixing forces in 
detail. 

To fix ideas we call B phase a "metal" and A phase 
the "vacuum" . Accordingly we put ua = n^A ~ ^ 
/a = 0. These last conditions do not change the solution 
but make the interpretation more transparent. 

Since in this case the number of particles in each phase 
is fixed (zero for the vacuum) we can minimize the energy 
per particle E = f /n given by: 



TIB 



i{x,nB) 



(Bl) 



This has to be minimized respect to the volume fraction 
taking into account that the density hb is also a function 
of the volume fraction given by the constrain hb = n/x. 
By putting the derivative respect to x of Eq. (Bl) equal 
to zero we obtain: 



Pb 



dx 



2 ^ni 

3 X 



(B2) 



The left hand side originates in the first term in Eq. (Bl) 
and is the intrinsic pressures of the metal, i.e the pressure 
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that the metal exerts on the surface. The right hand 
side is the pressure that the mixing forces, considered as 
"external" to the inhomogeneity, exert on the metal. We 
call this the mixing pressure (pm)- In equilibrium both 
pressures balance {pb — Pm ) • 

The mixing pre ssur e has two terms, the first [right 
hand side of Eq. (B2)] comes from the explicit depen- 
dence of the mixing energy on the volume fraction at 
constant ns and is proportional to u'{x). For an ide- 
ally symmetric u{x) (see Fig. ^ this term is positive for 
X < and is negative for a; > 0. We can say that this 
term tends to "compress" the metal in a less than half 
filled cell whereas it tends to "stretch" the metal (nega- 
tive pressure) in the opposite case. This is just the ex- 
pected tendency of the mixing energy to favor the closest 
uniform phase (cc = 1 or x = 0). The second term is due 
to the dependence of the mixing energy on the volume 
fraction through at constant particle number. An 
expansion of the cell at constant particle number pro- 
duces a decrease on which reduces the mixing energy 
[Eq. (|l^)]. This produces a negative pressure contribu- 
tion proportional to —u{x)/x. The net contribution is 
given by u' (x) - 2u{x) / {3x) [Eq. (|b|)]. It follows that for 
more that half filled cells the metal is subject to a net neg- 
ative pressure and for less than half filled cells the metal 
is subject to negative or positive pressures depending on 
the geometry and the volume fraction. The upper curves 
in Fig. H are proportional to u'{x) — 2u{x) / {3x) — 1. For 
drops the mixing pressure is positive for small x and then 
becomes negative whereas for layers the mixing pressure 
is negative for all x. 

The appearance of negative pressures indicates that 
the metal can be stable at densities which in the absence 
of LRC forces would be unstable so it is an indication 
of the stabilization effect of the LRC forces. In Fig. |lj 
we show E{nB) for a parabolic free energy. The intrinsic 
pressure (oc dE/duB) is negative for ns < and is 
positive for Ub > n'g. We will show below that stable 
solutions can be found in the region ub < n'g which are 
inaccessible (unstable) according to MC. 

The following example clarifies the physical mining 
of the negative pressures. Consider a neutral liquid 
with short range attractive forces. At negative pressure 
molecules will be at distances larger than the equilibrium 
distance and this implies an energy cost proportional to 
the volume. The system can relax by creating a surface 
and relaxing all molecules to the equilibrium distance. 
The energy cost proportional to the surface is much less 
than the energy gain proportional to the volume and this 
produce the MC instability. In the presence of meso- 
scopic frustrated PS this can not be done because for the 
drops the surface is not any more negligible respect to 
the volume. In fact the optimum drop ratio can be seen 
as the length scale at which successively breaking large 
drops subject to the negative mixing pressure is not any 
more convenient due to the surface energy cost. 

In the following we illustrate the behavior of the so- 
lution performing a graphical minimization of the en- 



ergy for the drop geometry and the parabolic free en- 
ergy Eq. ( |l8| ) . Instead of minimizing with respect to the 
volume fraction we use the constraint to eliminate the 
volume fraction in favor of ub {x — n/uB)- The energy 
per unit particle is given by: 
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(B3) 
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The /(f term has been eliminated with the MC condi- 
tion. The first term in the curly brackets is the bulk 
energy contribution. The mixing energy per particle is 

1/3 

em/(xnB) ~ u{x)/{xn^ ) and contributes to the last 
term in the curly brackets. The geometric factor u{x)/x 
gives the term in the square brackets. 

The equilibrium density is found by minimizing 
Eq. ( [B3D with respect to ub- In Fig. |l^ we show E — ijP 
as a function of ub for A — 0.3 and different values of 
n' . The thick line is the energy of the uniform metal [the 
first term in the brackets in Eq. ( p33| )] and is minimized 
at the MC density Uq . 
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FIG. 14. Normalized energy per particle as a function of 
TiB/n% for A — 0.3. The thick line correspond the uniform 
phase and the thin lines to the drop state with n' changing 
from zero (top) to one (bottom) in steps of 0.1. The crosses 
indicate the drop solution and the uniform solution. 

For very small n (or x) the geometric factor is constant 
and the mixing energy contribution goes as I/hb^^^- 
This shifts the minimum to values of the density lar ger 
than Uq as can be seen from the upper curves of Fig. 04 
where the energy per particle is given for various values 
of the global density and A = 0.3. This is due to the 
positive pressure exerted by the mixing energy of drops 
at small volume fraction and explain the behavior of the 
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ub density in the limit n' — > (Fig. |^). As the density 
increases the density dependence of the geometric factor 
tends to reduce the minimum to lower densities. 

As a by product this computation illustrates the sta- 
bilization of a uniform solution by the long range inter- 
action and the first order like nature of the transition. 
Above n' ~ 0.6 the uniform solution becomes suddenly 
more favorable (see also Fig. ^ . Notice that this density 
is well inside the MC coexistence region (0 < n' < 1) 
showing the uniform solution stabilization effect. 

It is important to remark that the whole behavior can 
change if the surface energy a had a strong density de- 
pendence. For this reason the interpretation of B phase 
as a metal should be taken with caution since in general 
in a metal the surface energy will depend strongly on 
density. Specific examples will be treated in paper II. 
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